Continuum percolation of simple fluids: Energetic 
connectivity criteria. 

Luis A. Pugnaloni 
Procter Department of Food Science, University of Leeds, 
Leeds LS2 9JT, UK 

Ileana F. Marquez and Fernando Vericat 
Instituto de Fisica de Liquidos y Sistemas Biologicos 
(IFLYSIB)-UNLP-CONICET cc. 565 - (1900) La Plata, Argentina 

February 1, 2008 

Abstract 

During the last few years, a number of works in computer simulation 
have focused on the clustering and percolation properties of simple fluids 
based in an energetic connectivity criterion proposed long ago by T.L. Hill 
[J. Chem. Phys. 23, 617 (1955)]. This connectivity criterion appears to 
be the most appropriate in the study of gas-liquid phase transition. So 
far, integral equation theories have relayed on a velocity-averaged version 
of this criterion. We show, by using molecular dynamics simulations, that 
this average strongly overestimates percolation densities in the Lennard- 
Jones fluid making unreliable any prediction based on it. Additionally, 
we use a recently developed integral equation theory [Phys. Rev. E 61, 
R6067 (2000)] to show how this velocity-average can be overcome. 

1 Introduction 

Clustering and percolation in continuum systems are concepts of great interest 
in chemical physics. Phenomena such as nucleation |l], hydrogen bonding 
conductor-insulator transitions sol-gel transitions and bridging in gran- 
ular materials are frequently studied in terms of clusters and percolation 

i- 

The first application of statistical mechanics formalism to describe clustering 
in equilibrium classical systems was done by Hill . In Hill's theory, the concept 
of cluster is directly related to the idea of bonded-pairs. A bonded-pair is a set 
of two particles that are linked by some direct mechanism. Then, a cluster is a 
set of particles such that any pair of particles in the set is connected through 
a path of bonded-pairs. We call this clusters "chemical clusters" to distinguish 
them from the non-pair-bonded clusters we have introduced in a previous work 
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[||. A system is in a percolated state if it contains a cluster that spans the 
system. 

From Hill's theory, we see that, in order to decide when two particles are 
bonded, a connectivity criterion is necessary. This connectivity criterion has 
to be defined according to the phenomenon under study ||^. In the search for 
clusters of atoms in a monatomic gas that marks the onset of a phase transition, 
Hill proposed a simple energetic criterion: two particles are bonded if their 
relative kinetic energy is less than minus their relative potential energy Q. 
However, this criterion is difficult to implement from a theoretical point of view 
and simpler ones were preferred. The very first simplification was set by Hill 
himself. Instead of using his general criterion, for which the relative positions 
and velocities of the relevant pair of particles come into account, he averaged 
out all the possible velocities of the particles. Then, a velocity-averaged (VA) 
criterion can be stated as follow: two particles are bonded with probability P{r) 
if their relative potential energy is negative. P{r) is calculated as an average 
probability over all possible relative velocities of a pair of particles . 

Coniglio et. al. have obtained an Ornstein-Zernike-type relationship for the 
pair connectedness function g'l'(ri, r2). This function is proportional to the joint 
probability density of finding two particles belonging to the same cluster and at 
positions ri and Y2, respectively Then, by integrating g^(ri,r2), the mean 
cluster size S and the percolation density pp (i.e. the value of p for which S{p) 
diverges) can be obtained. 

It is worth mentioning that most of the theoretical studies on connectivity 
and percolation in continuum systems (see for example Ref. jTl], p^ , p^ ) were 
focused in the rather simple Stillinger's connectivity criterion | ]14| . It states 
that two particles are bonded if they are separated by a distance shorter than 
a given connectivity distance d. In this case, d is an ad hoc parameter, which 
must be chosen on physical grounds. Although this criterion might be sensible 
in the study of certain insulator-conductor transitions, it is unrealistic regarding 
clustering in saturated vapours. 

More recently. Hill's criterion has been reconsidered in molecular dynam- 
ics studies of small clusters ^ and the critical percolation behaviour of 
Lennard- Jones fluids It has been suggested that the percolation line — 

the line that separates the temperature-density phase space into percolated and 
non-percolated states — might be experimentally observable More- 
over, cluster analysis based in Hill's criterion seems to be useful in locating 
gas— liquid coexistence curves [0. Notice that molecular dynamics simulations 
are mandatory if Hill's criterion is used to identify clusters since Monte Carlo 
algorithms do not provide the velocities of the particles to decide whether the 
relative kinetic energy of the a given pair outweighs its relative potential energy. 

Since Coniglio's theory deals only with the positions of the particles. Hill's 
criterion cannot be achieved. Instead, the VA criterion was used by Coniglio et. 
al. [ p^ to analytically calculate, for a potential made up of a hard core plus an 
attractive interaction, the percolation loci in a crude mean-field approximation. 
Although they obtained that the percolation line ended just below the critical 
point, we show that this is an artefact of an additional approximation they 
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introduced, which is not very suitable for the low temperatures involved. 

In this work we present the results of a molecular dynamics (MD) inves- 
tigation of the percolation loci for a Lennard- Jones system by using the two 
previous energetic criteria, i.e. VA and Hill's energetic criteria. We add a max- 
imum connectivity distance d to the criteria to avoid unrealistic bonding. Two 
particles separated by a distance greater than this maximum distance are con- 
sidered as non-bonded even if their relative velocity is by chance zero . We 
show that the velocity average introduces an important overestimation on the 
predicted percolation density at all temperatures. Therefore, the approximation 
introduced by this velocity average is not sensible. As a consequence, a way to 
implement the full Hill's connectivity criterion in an integral equation approach 
becomes essential for the future of these theories. Here, we show how a recently 
developed integral equation theory for continuum percolation [20[| can be used 
to solve this problem. 

The rest of the paper is organised as follow. In Sec. 2 we describe the MD 
simulations and show the unsuitability of the velocity average in the energetic 
bonding criterion. In Sec. 3 we briefly discuss Coniglio's percolation theory 
and two examples of its application to the VA criterion. In Sec. 4 a generalized 
integral equation theory is used to illustrate how the full Hill's energetic criterion 
can be fitted into this formalism. We conclude with some consideration about 
the possible future developments of this approach. 



2 Molecular dynamics 



Molecular dynamics simulation of 500 particles in a cubic box with periodic 
boundary conditions in the NVT ensemble was performed by a leap-frog al- 
gorithm with velocity correction ||2l| . Particles interact via the Lennard- Jones 
potential 

'(T\12 /(T\6l 



v{r) ~ 4e 



(1) 



Time step was chosen as At* — IS.ta~^ ^ksT j (em) = 0.01. Quantities are 
averaged over lO'^ configurations picked up each 50 At after stabilization. A cut 
off distance equal to 2.7(t was used in the pair potential. 

We identify the clusters by using both, the Hill's criterion for the bonded 
probability defined as 0: 

J 1 pf j/4m < ~v{vi,v^) andri,j < d 
^^^''i'P''3> - \ p2^./4m ^ ~v{r„ r^) or > d ^ ' 

(jij and Pi j are the relative position and momentum of particles i and j), and 
the VA criterion, which correspond to |Q, |l^ : 

p , x_ / w(r,j) > Oorr.j > d 

^av{r,,j) - I ^j3/2^ -t;(r,,,)/fcBT]/r[3/2] «(r,,,) < Oanrfr,,, < d 

(T[a] is the gamma function and 7[a,x] is the incomplete gamma function). 
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In the last case, for each pair of particles that satisfy v{rij) < and r,;.j < d, 
we generate a randomnumber z, between and 1. If z < ■^[3/2,—v{rij)/kBT]/r[3/2] 
we consider that the particles form a bonded-pair, otherwise we do not. Note 
that this criterion can also be used in Monte Carlo simulations because it does 
not need the knowledge of the particle velocities. We also identify clusters 
according to Stillinger's criterion for the sake of comparison. To identify the 
clusters from the list of bonded-pairs we use the Stoddard algorithm ^ . 

A system is said to be in a percolated state if a cluster that spans through 
the replicas is present 50 percent of the time Then, a percolation transi- 
tion curve, which separates the percolated and the non-percolated states of the 
system, can be drawn above the coexistence curve in the T — p phase diagram. 

In Fig. 1 percolation loci for Hill's, VA, and Stillinger's connectivity criterion 
are presented for d — 2a. The coexistence curve obtained by Panagiotopoulos 
[ p^ is also included for comparison. Theoretical predictions from connectedness 
integral equation theory are shown for VA and Stillinger's criteria only (see 
Sec. 3). The two distinctive differences of the results for the two energetic 
criteria, as compared with the Stillinger's criterion, are: i) a strong dependence 
on temperature of the percolation curves and ii) a shift towards higher densities. 

As it can be seen in Fig. 1, the VA criterion yields higher percolation 
densities than Hill's criterion. As it is expected, the approximation due to 
the average is less important as temperature drops: for low temperatures, the 
kinetic energy is small compared to the potential energy and an average over 
velocities is a good approximation. Nevertheless, the overestimation of the 
percolation density is never small. Therefore, it is expected that any theory 
based on the VA criterion will fail in describing the proper percolation (and in 
general clustering) behaviour of simple fluids. Any success should be attributed 
to cancellation of errors rather than accuracy of the theory. Indeed, we shall see 
in the next section that a Percus-Yevick-like approximation in the connectedness 
Ornstein— Zernike integral equation yields, for the VA criterion, an apparently 
good percolation line as compared with the full Hill's criterion. 

It is worth mentioning that the results for the Hill's criterion are similar to 
the ones obtained by Campi et. al. in a more extensive MD simulation. 
However, they obtained a faster increase in the percolation density with increas- 
ing temperature. The differences between our results and those from Campi et. 
al. could be related to the small size of our system and most importantly to 
the criterion for deciding whether the system is percolated or not. In Ref. |0 
the authors considered that the system is on the percolation line if the second 
moment of the cluster size distribution n'{s) that excludes the largest cluster 
reaches its maximum. Nevertheless, we are interested in the suitability of the VA 
bonding criterion and we expect that any criterion to calculate the percolation 
line might be affected to a similar extent by this simplification. 
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3 Continuum clustering theory for non- velocity- 
dependent bonding criteria 

3.1 General Theory 

In order to study clustering in a system composed of N classical particles in- 
teracting via a pair potential t'(ri,r2), Hill separated the Boltzmann factor 
e(ri,r2) = exp[— /3w(ri, r2)], into bonded (f) and unbounded (*) terms [Q: 
e(ri,r2) = et(ri,r2) + e*(ri,r2). As usual f3 — l/fc^T, being fc^ the Boltz- 
mann constant. Since e'l'(ri,r2) represents the basic probability density that 
two particles at position ri and r2 are bonded, this separation yields a diagram- 
matic expansion for the partition function in terms of "chemical" clusters. For 
convenience we express Hill's separation as follow: 

e'l'(ri,r2) = P(ri, r2) exp[-/3'(;(ri, r2)] (4) 



e*(ri,r2) = [1 - P(ri, rj)] exp[-/3i;(ri, r2)] (5) 

where P(ri,r2) is given by Eq. (H) in the case of the VA energetic criterion. 

Fugacity and density expansions have been found, within Hill's formalism, 
by Coniglio and co-workers for the pair connectedness function g''{ri,r2). 
This function is proportional to the joint probability density of finding two 
particles belonging to the same cluster and at positions ri and r2, respectively. 
Moreover, by collecting nodal and non-nodal diagrams in these expansions an 
Ornstein-Zernike-type relationship is obtained 

5^(i"i,i"2) = c^(ri,r2) +p y"c^(ri,r3)g^(r3,r2)dr3. (6) 

The function c^(ri,r2) is the direct pair connectedness function. By posing 
a closure relation between this function and g^(ri,r2), Eq. (||) can be solved. 
Then, the mean cluster size can be obtained as 

Sip) = 1 + j^^^ j g\v,,V2)dv^dv2 (7) 

were the integration is carried out over the whole system volume. 

Percolation of an ensemble of particles is concerned with the existence of 
clusters that become macroscopic in size. Therefore, the percolation transition 
occurs at a critical percolation density pc which mathematically verifies 

liW^S{p) = oo. (8) 

In the rest of this section we present results for two systems were the VA 
criterion is applied. 
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3.2 A simple example 

A theoretical calculation of the percolation loci for the VA criterion [see Eq. 
(||)] was done by Coniglio et. al. for a system in which the pair potential 
has the following form: 

f oo ri,2 < To 

i;(ri,r2) = w(ri,2) - <{ ^. f ro V „ ^„ • (9) 



However, they did not use a maximum connectivity distance d to avoid unphys- 
ical connections. Moreover, they used an expansion for the incomplete gamma 
function up to first order. We have recalculated numerically the percolation 
curve using d — 2a and the complete expansion of the gamma function. The 
corresponding expression for the percolation density pp is fl^ : 

where c^(fc) is the Fourier transform of the direct pair connectedness function. 
We have used the closure relation for Eq. proposed by Coniglio et. at, i.e. 

ct(ri,r2) = /t(ri,2) (11) 

where /''{rij) = e^rij) = exp[— f}v{rij)]Pav{rij) is the basic probability for 
finding two particles separated by a distance rij and bonded. The function 
f^rij) is called the bound Mayer function. 

In Fig. 2 we reproduce the coexistence curve and the percolation line, for 
the VA criterion, from Ref. The corrected results, using the complete 

expansion for the incomplete gamma function and the maximum connectivity 
distance d, are also included. Additionally, we show the percolation line for 
the full Hill's criterion (see Sec. 4). As we can see, using the full expansion 
of the incomplete gamma function is important since results are very much 
affected. The importance of the maximum connectivity distance d is negligible. 
We calculated the percolation line for increasing values of d — up to infinity — 
without obtaining significant changes. 

Based on the earlier results, it was suggested Q that this theory is rather 
satisfactory considering its simplicity because the predicted percolation line 



ended close to the critical point as MD simulations show 1 17 for the full Hill con- 
nectivity criterion. However, this might be due to a cancellation of errors. On 
the one hand, the actual percolation line for the VA criterion should be shifted 
to even higher densities as it can be seen in Fig. 1 for the Lennard- Jones sys- 
tem. On the other hand, the approximation for the incomplete gamma function 
seems to affect the result largely. 

3.3 Numerical results for the Lenard-Jones system 

In order to solve Eq. ^ for the Lennard-Jones potential we have implemented 
Labik's numerical algorithm pa. We have already seen that the simple closure 
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relation given by Eq. (|11|) is very crude. Hence, we use the more reliable closure 
available, i.e. the Percus-Yevick-like relation 

gHna) - [/*(^i,2) + l][.9^(r-i.2) - c\n^2)] + exp[f3v{n^2Mri,2)fHri,2). (12) 

In Eq. /*(ri,2) = e*(ri,2) - 1 = exp[-/3z;(r,j)][l - Pav{n^j)] - 1 is the 

unbound Mayer function. We recall that Pavifi.j) is given by Eq. (|) for the 
VA energetic criterion. 

By using Labik's algorithm we solve the coupled system of equations (^ 
and (|l2|). Then, mean cluster size can be obtained for given T and p. For a 
given temperature, percolation density is determined by fitting the power law 
5*^^ oc \p ~ pc\'^ for p ^ Pc p3| . Although critical exponents obtained by fitting 
iS* as a function of p are difficult to obtain accurately in these cases, percolation 
densities are very reliable. 

In Eq. ( |T^ ) the "thermal" pair correlation function g{ri,2) is involved. We 
have obtained 5(^1^2) by numerically solving the "thermal" Ornstein-Zernike 
equation closed with the "thermal" Percus-Yevick approximation pof . 

In Fig. 1 we show the results for the VA criterion and the Stillinger's cri- 
terion for d = 2(7. As we can see, the results from Percus-Yevick theory agree 
rather well with MD simulations for Stillinger's criterion. However, for the VA 
criterion, the theory is unsatisfactory. Nevertheless, we see a good advance with 
respect to the very low percolation densities predicted in the previous subsec- 
tion, for a similar system, by the closure relation (pi]). We believe that the main 
cause of the failure for the VA criterion is related to the "input" function g{ri^2)- 
The value of this function is less well predicted by the "thermal" Percus-Yevick 
approximation at these "high" densities and temperatures. Notice that these 
types of theories were developed to describe the liquid state, particularly near 
the triple point. At low densities, as the involved for the Stillinger's criterion, 
these theories are valid even for temperatures above the critical point. 

Although the previous theory seems to be inappropriate to describe the 
percolation transition for the VA criterion, it seems to agree better with the 
MD result for the full Hill's criterion. This is, of course, due to a cancellation of 
errors between the velocity average and the Percus-Yevick approximation. Since 
Hill's criterion produces a percolation line at intermediate densities, we should 
be able to use the Percus-Yevick approximation successfully to study this case. 
However, getting rid of the velocity average is the main task to accomplish. In 
the following section we show how this can be achieved. 

4 Continuum clustering theory for velocity-dependent 
bonding criteria 

4.1 General theory 

In this section we summarize the main results of a theory we have presented 
elsewhere M] to describe the clustering and percolation for bonding criteria 
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that involve the velocity of the particles as well as their positions. 

For a system of N classical particles interacting via a pair potential v{ri, rj) 
we define a density correlation function p(ri, r2, Pi, P2) which is N{N — 1) times 
the probability density of finding two particles at phase space points (ri, pi) 
and (r2, P2) respectively: 



P(ri>r.,Pl,P2) = h^N^^^Q^^y^T^ (13) 

X / Y[exp[-(3^]l[l[eM-Mn,r,)]dr''''dp^'\ 

i=l i=lj>i 

Here h is the Planck constant and Qn{V,T) the canonical partition function 
of the system. Then, in the same spirit of Hill and Coniglio et. al. |l^ , we 
separate exp[— /3w(ri, rj)] into connecting and blocking parts 

exp[-/3w(r,,rj)] = /t(r,, r^, p,, p^) + /*(ri, r^, p^, p^) + 1 (14) 

Here /'''(r^, r^, p^, Pj) represent the basic probability density that two particles 

at configuration (r^ , rj ,Pi,Pj), are bonded. The shorthand notation /''' (r^ , rj , p,j ,pj) = 

fij ('^ ~ ^^^^ sometimes used. 

Substitution of Eq. (|l|) in Eq. (|l|) yields 



p(ri,r2,Pi,P2) = ;,3iVjv!Q^(y,T) ''''P^^''^''^''''^^ ^^^^ 

/■ ^ 2 

1=1 

where the sum is carried out over all possible arrangements of products of func- 
tions and fli- 

It should be noted that the functions //^ and f*j can depend on momenta 

as well as on the positions of the two particles, but the sum of f^^ and f*^ must 
be momenta independent in order that Eq. (^4|) be satisfied. Except by this last 
condition, the functions fjj and f*j are otherwise arbitrary for thermodynamic 
purposes. Obviously, we choose them in such a way that the Hill's definition of 
bonded particles is achieved: 

flj = exp[-/3^^(r,j)]P(r,j, p,,j) (16) 



f*^ = exph/3i;(r,,,)][l - P(r,,„p,,,)] - 1 (17) 

where P(rij, p^j) is given in Eq. (||). 

Each term in the integrand of Eq. ( |l5| ) can be represented as a diagram 
consisting of two white ei and e2-points, N ~ 2 black e^-points and some 
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and f*j connections except between the white points. Here we take = 

2 

exp[— /3|^]. White points are not integrated over whereas black points are 
integrated over their positions and momenta. AU the machinery normaUy used 
to handle standard diagrams in classical liquid theory can be now extended 
to treat these new kind of diagrams. By following Coniglio's recipe to sep- 
arate connecting and blocking parts in correlation functions i.e. (7(ri,r2) = 
5i'''(ri, r2, Pi, P2) + 5*(ri, r2, Pi, P2), we obtain an Ornstein-Zernike like integral 
equation for .g'f(ri, r2, pi, P2) 

5rt(ri,r2,pi,P2) = ct(ri,r2,pi,P2) (18) 
+ y P(i"3, P3)c'''(ri, ra, Pi, P3)5'^(r3, r2, P3, P2)rfr3dp3. 

Here p(ri, pi)p(r2, P2)g^(ri, r2, pi, P2) is iV(iV — 1) times the joint probabil- 
ity density of finding two particles at positions ri and r2 with momenta pi and 
P2, respectively, and belonging to the same cluster, where the bonding criterion 
is given by Eqs. (^6|) and ([l7|); being 



p(ri,Pi) = jy ^ ^ y p(ri,r2,pi,P2)dr2dp2. (19) 

The function c^(ri, r2, pi, P2) denotes the sum of all the non-nodal diagrams 
in the diagrammatic expansion of (7^(ri, r2, pi, P2). We remember that a nodal 
diagram contains at least a black point through which all paths between the two 
white points pass. 

For homogeneous systems we have 



5'^(i-i,i-2,Pi,P2) = ct(ri,r2,pi,P2) + TT^f^jH (20) 



J exp[-(3^]c'{ri,r3, pi, P3).g^(r3, r2, ps, P2)dr30?p3 



To get an integral equation from Eq. (|l^) is necessary a closure relation 
between g''{ri, r2, pi, P2) and c^^(ri, r2, pi, P2). As an example, we can take the 
Percus-Yevick approximation ^(ri, r2) exp[f3v{ri, r2)] = 1 + A^(ri, r2), where the 
function iV(ri, r2) is the sum of the nodal diagrams in the expansion of g{ri, r2). 
Separation into connecting and blocking parts i.e. (7(ri, r2) = 5^(ri, r2, pi, P2) + 
g*(ri,r2,pi,P2) and 7V(ri,r2) = 7Vt(ri, r2, pi, P2) + iV*(ri, r2, pi, P2), yields 

ff'^(ri,r2,pi,P2) = [/*(ri,r2,pi,P2) + 1] [g^(ri, r2, pi, P2) - c^(ri, r2, pi, P2)] 
-I- exp[/3w(ri, r2)]5(ri, r2)/'l'(ri, r2, pi, P2). (21) 



Eq. ( p^ closed by Eq. ( [2l| ) gives an integral equation for g^(ri, r2, pi, P2). 
From the function 5^(ri, r2, Pi, P2) we define the pair connectedness function 

5^(i"i,i"2)= y"p(ri,pi)p(r2,p2).g^(ri,r2,pi,p2)dpidp2. (22) 
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This function is the joint probabiHty density of finding two particles that be- 
long to the same cluster, within the Hill's criterion, at positions ri and r2, 
respectively. 

Then, the mean cluster size S is given by 



S = 1 



(N- 



^Tyy j g\vi,T2)dvidv2. 



(23) 



4.2 A simple example 



By using the previous extension of Coniglio's work, we have calculated the 
percolation loci for the full Hill's criterion. We have used the pair potential of 
Eq. . In this case we solve Eq. with a closure relation equivalent to Eq. 
TTl), i.e. 



(24) 



c^(ri,r2,pi,P2) = /^(ri,r2,pi,P2). 
where /'''(ri, r2, pi, P2) is given by Eq. (|l6|). Then 







ci'(ri,r2,pi,P2) = 



exp 



ri,2 < ro or {piaY /Am > -w(ri,2) 
ri,2 > ro and {pi^2)^ /4:m< ~v{ri^2) 

d25) 
) as 

follow. Firstly, taken into account the homogeneity of the system, we perform 
a Fourier transform on Eq. (|20|) with respect to ri^2 



5'''(ki,2,Pi,P2) = c'l'(ki,2,Pi,P2 



(27rmfci3r)3/2 



(26) 



where 



X J exp[-/3^]c1'(ki,2,Pi,P3)ff^(ki,2,P3,P2)dp3 

g'^(ki,2,Pi,P2) = y g'^(ri,2,Pi,P2)exp[-iri,2ki,2]dri_2 (27) 
and an analogous expression for c^(ki^2, Pi, P2)- 

2 

By evaluating the previous equation in ki.2 = 0, multiplying by exp[— 
and integrating with respect to P2 we obtain 



glipi) ^ 4(Pi) + 



with 



f P^ 

/ exp[-/3^]4(pi,p3)gj(p3)dp3 
ffJ(Pi,P2) = g'''(ki,2 = 0,Pi,P2), 



(28) 



(29) 
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r p2 

5o(Pi) = / exp[-/3^]5j(pi,p2)cip2, (30) 

and analogous expressions for ct(ki_2, Pi, P2)- 

Because cj 2 in Eq. (|2^) depends just on the modulus of ri — r2 and Pi — P2, 
we can perform a Fourier transform with respect to pi yielding 

9li^) - [-^ ) exp[-^]ct (.) + (,^J^^)3/2 4(^)g^(^) (31) 
where we have defined 

= J gi(pi)exp[-iwpi]dpi (32) 

f p^ 

9o{'^) ^ / ffo(P3)exp[-/3^]exp[-icjp3]dp3 (33) 

= J 4)i\Pi - P3|)exp[-iwpi_3]dpi,3. (34) 
From Eq. (^), the mean cluster size can be expressed as 

5* = 1 + \ y"p(ri,pi)p(r2,P2)g^(ri,r2,pi,P2)dpidp2drirfr2, (35) 
or, as a function of gl{uj), 

+ I ^,3 go(^ = 0)- (36) 
[znmKBi ) 

Taking .gj(a;) = f ^X^) 6xp[— ^^^^]cJ((jj) as a first order approximation 
and iterating Eq. ( |3l|) up to convergence we can calculate the mean cluster size 
for given p and T. 

As in the previous section, percolation density is determined by fitting the 
power law S^^ cx jp — f^'' P ^ Pc- 

In Fig. 2, we present the percolation line obtained for Hill's criterion. As 
for the VA criterion, the simple closure relation implemented yields unrealistic 
percolation densities. However, it can be seen that the percolation line is shifted 
to the left, as compared with the one corresponding to the VA criterion. This 
line is also less temperature dependent — it presents a sharper slope in the T — p 
plane. These trends are consistent with the MD simulation for the Lennard- 
Jones fluid (see Fig. 1). The discrepancy between both criteria is negligible 
only at low temperature as expected. 

Although a satisfactory theory for the continuum percolation of a simple 
fluid within the Hill's energetic connectivity criterion is still due, we believe the 
main step forward has been presented here. The next attempt to improve the 
theory should be the implementation of the Percus-Yevick connectedness closure 
relation (pl|) . This can be done by solving numerically the coupled system (EOh 
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and This is a non-straightforward task: tensors of rank 9 must be used 

to storage correlation functions and a double convolution has to be solved. 
Recently, F. Lado |^ has developed a method to calculate the pair correlations 
for a polarizable Lennard- Jones system which is based in an expansion of the 
relevant functions in orthogonal polynomials. The integral equation he has to 
solve is similar to Eq. (|l^) but with pi and P2 representing the induced dipolar 
moments of the particles 1 and 2 respectively. A modification of this method 
might be used for the clustering problem of Lennard-Jones particles by using 
the Percus-Yevick like closure of Eq. . 



5 Concluding remarks 

In this work we have shown that the common introduction of an average over 
velocities in the connectivity energetic criteria for simple fluids strongly overes- 
timates the percolation density. Moreover, we have shown how this average can 
be avoided by using a more general continuum percolation theory. Although 
this theory was initially developed to implement time-dependent connectivity 
criteria, it is specially suited to the Hill's energetic criteria. 

Even though we have applied the new formalism to avoid the velocity- 
average, we used a very simple approximation to close the generalized Ornstein- 
Zernike connectedness equation. Nevertheless, the results are promising in that 
the trends obtained by MD simulation are reproduced. A more sophisticated 
closure relation should be implemented in order to obtain a more satisfactory 
theory. 
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Figure Captions 

Figure 1: Coexistence and percolation curves for the Lennard- Jones fluid: 
MD (symbols) and Percus-Yevik connectedness theory (lines). Densities and 
temperatures are in units of and e/ks, respectively. Solid diamonds (with 
trend line) correspond to the coexistence curve. Percolation loci for Stillinger's 
criterion (open squares and dotted line). Hill's criterion (open circles), and VA 
criterion (open triangles and solid line) are presented. 

Figure 2: Coexistence and percolation lines for a simple fluid. Full and 
dashed lines correspond, respectively, to the percolation loci for Hill's criterion 
[calculated from Eqs. (|o|) and (|2^)] and VA criterion [calculated from Eqs. 
(||) and (|ll|)]. Dotted hue for VA criterion is reproduced from Ref. ||l^ to be 
compared with the "corrected" results (dashed line). The coexistence line is 
also reproduced from Ref. |19|] . The particles of the fluid interact through the 
pair potential given by Eq. (O) . Densities and temperatures are reduced with 
the corresponding critical values. 



14 



